Periodic Orbits in a Simple Ray-Splitting System 



Debabrata Biswas * t 
Center for Chaos and Turbulence Studies, Niels Bohr Institute 
Blegdamsvej 17, Copenhagen 0, Denmark 



Abstract 



We study ray dynamics in a square billiard that allows mode conversion 
and is parametrised by k, the ratio of the velocities of the two modes. At 
k — > 1 + , conversion occurs at every reflection and periodic orbits proliferate 
exponentially. As k increases beyond y/2, the collection of daughter rays 
explore only three momentum directions and mode conversion is progressively 
inhibited. We provide an algorithm for determining periodic orbits when 
K > v2 and show numerically that exponential proliferation persists around 
\/2 but as k increases, a crossover to sub-exponential behaviour occurs 
for short periods. We discuss these results in the light of conservation laws. 
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I. INTRODUCTION 



Billiards are useful examples of dynamical systems and display the wide variety of phe- 
nomena associated with Hamiltonian flows. The ray equations commonly considered follow 
from a short wavelength expansion of the Schrodinger equation with Dirichlet (or Neumann) 
boundary conditions. A particle thus moves freely between collisions at the boundary where 
it suffers specular reflection and one can observe regular or chaotic motion depending on the 
shape of the boundary. As an example, a square billiard generates regular dynamics that is 
restricted to a torus in phase space due to the existence of two well behaved constants of 
motion. In contrast, generic trajectories in the stadium shaped billiard explore the entire 
constant energy surface and hence the system is ergodic. Moreover, orbits that are nearby 
initially move exponentially apart with time and hence the system is said to be metrically 
chaotic. 

These differences also show up in the proliferation rates of periodic solutions. In case 
of regular motion, periodic orbits exist in 1-parameter families and their number increases 
quadratically with time period. In contrast, chaotic dynamics is accompanied by an expo- 
nential proliferation of periodic orbits, a phenomenon referred to as topological chaos. 

The manner in which periodic orbits organise themselves in closed systems is strongly 
linked to the existence of sum rules arising from conservation laws. For example, the fact 
that a particle never escapes implies that for chaotic systems 



where the summation over p refers to all primitive periodic orbits, T p is time period, J p is the 
stability matrix evaluated on the orbit and the symbol (.) denotes the average value of the 
expression on the left. Since the periodic orbits are unstable and isolated, |det(l — J/)| — 
e AprTp , where X p is the Lyapunov exponent of the orbit. The exponential proliferation of 
orbits is thus implicit in eq. (jl|). 

For the square billiard though, the appropriate sum rule is |2|.|3|] : 
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where v is the velocity, Tm,n is the time period of a periodic orbit with winding numbers 
(M,N) and A is the area of the billiard. The quadratic law for the number of periodic orbits 



In both cases, the stability of periodic orbits leads to the respective proliferation laws. 
However, there can exist situations where rays split up at the boundaries of billiards and 
thus suffer a decay in intensity. The basic conservation law then demands that periodic 
orbits should proliferate faster than usual to compensate this loss. 

We explore here the case of a square billiard that admits ray-splitting and is parametrised 
by k, the ratio of velocities of the two modes. For n — > 1 + , mode conversion occurs at every 
reflection from the boundary and daughter rays multiply as 2 k where k denotes the number 
of reflections. Moreover, the collection of daughter rays (from a single parent) can access an 
increasing number of momenta directions with time. As k increases, the range of angle in 
which conversion can occur at adjacent edges decreases and consequently the average number 
of daughter rays produced with every reflection decreases. For k > \/2, the possibility of 
mode conversion is restricted further and daughter rays can access only three momentum 
directions. Simultaneously, the range of angles at which mode conversion can occur shrinks 
as k increases beyond n c = In the limit k — > oo conversions are not allowed and the 
system is then a normal billiard without ray-splitting. 

In the setting described above, we study the proliferation law of periodic paths. In 
particular, we explore the region k > \/2 and provide an algorithm for determining periodic 
paths. We find that for k = 1.429 (> k c ), the proliferation law remains exponential but as 
k increases further, mode conversion is inhibited and the proliferation law shows a crossover 
to sub-exponential behaviour for short periods. 

Before discussing these in detail, it is important to note that periodic orbits form the 
skeleton on which modern semiclassics is built [|j] and are used to understand the spectrum 
in quantum systems. In elastodynamics where mode conversion does occur, studies of the 
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statistical properties of the resonance spectrum are influenced to a large extent by these 
considerations [§,§J. In fact, a fourier transform of the resonance spectrum provides evidence 
of the role of periodic orbits in elastodynamics |J. The present study is thus significant in 
that it provides the first systematic computation of periodic orbits and underscores its 
complex organisation even in simple geometries. 

In the following section, we describe in more detail the system that we study and analyse 



the dynamics of the converted rays. Section [TTj deals with the algorithm for determining 



periodic orbits and our numerical results are discussed in section [TV 



II. A SIMPLE RAY-SPLITTING SYSTEM 

Ray-splitting is a common phenomenon in geometrical optics and can be observed in 
several other situations |7]-|9[]. The propagation of elasto-mechanical waves in solids is an 
example that is of much current interest |5]|| and we study here the short wave-length limit 
of this problem for the isotropic case in two dimensions. For small displacements, the wave 
motion is governed by the Navier equation [1(J : 



pV 2 u + (A + p) VV • u = pd 2 u/dt 2 (3) 



where u is the displacement, A and p are the two Lame constants [JTTJ] and p is the density. 
It is common to express the displacement as a sum of two parts generated respectively by a 
scalar and a vector potential (u = u x + u 2 where u x = V0 and u 2 = V x B) so that eq. ([3]) 
separates into two second-order equations [flOfl : 



V 2 = \d 2 <\>jdt 2 (4) 
c P 

V 2 B = \d 2 B/dt 2 (5) 

c s 

where cp = (A + 2p)/p and c| = p/p. The medium thus has two natural velocities, cp and 
cs where the subscripts P and 5* refer to the pressure (longitudinal) and shear (transverse) 
waves respectively. The two waves interact only at the boundaries where they may suffer 
mode conversion. 



Consider for example an S or a P wave incident on a planar stress-free (tractionless) 
boundary. The reflected part consists in general of both an S and a P wave and the 



reflection law involving an S and a P wave is given by Snell's relation [Kj: 

cos(6 l p) = kcos(9s) (6) 

where the angles are measured with respect to the tangent at the collision point on the 
boundary and k = cp/cs- When the reflected and incident wave belong to the same type 
(both S or both P), the angle of reflection equals the angle of incidence. Note that eq. |6| 
implies the existence of a critical angle, 9c = cos _1 (l/ft). For 9s < 9c, no conversion can 
occur and the wave suffers only a specular reflection. 

The incident and reflected amplitudes are also determined by the boundary conditions. 
For S — > S + P process, and for the case where the S wave polarisation lies in the plane of 
incidence, the intensities carried away by the reflected S and P wave are fL0|: 
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sin(2fl s ) sin(2fl P ) + k 2 cos 2 (2# 5 ) v ' 

and Isp = 1 — lss respectively. Similarly, for the P — > S + P process, Ipp = lss an d 
Isp = I ps provided again that the S wave polarisation lies in the plane of incidence. In 
general, the full 3-dimensional problem requires a decomposition of the polarisation vector 
into components normal and parallel to the plane of incidence. The component in the plane 
undergoes the process described above while the normal component suffers no conversion. 

The small-wavelength limit of this wave motion restricted to a plane can be treated 
as the relevant ray-tracing problem in 2-dimensions. The geometry we consider here is 
a square billiard of length, L = n/2. In a situation where no mode conversion occurs, 
the square billiard is integrable and its periodic orbits increase quadratically with length. 
The corresponding Schrodinger equation can be solved easily and the quantum spectrum 
is described exactly by periodic orbits. The situation in elasdo-dynamics is however quite 
the opposite and the only resonances known analytically form a small fraction of the total 
number [ITJ. To the best of our knowledge, the ray-tracing problem in this geometry has 



not been studied before and we proceed to understand that now. 



The parameter range that is physically accessible is n > 1 and as an example we study 
the case k = 1.429. Fig. 1 illustrates the conversions that can occur in this system. For 
the sake of visualisation, it is convenient to consider an unfolded trajectory in the full plane 
generated by reflections of the fundamental domain about its sides (see Fig. 2) . A reflected 
trajectory without conversion thus continues through the boundary without any change in 
angle while a ray that has suffered conversion is equivalent to a refracted ray that undergoes 
a change in angle. A cell identical to the fundamental domain can then be labelled by the 
winding numbers (M, N) where 2ML is the displacement of any point in the cell along the 
X-axis and 2A^L along the Y-axis. Fig. 2 shows the cells with winding numbers (1,1) and 
(2,2). 

Consider then a family of S rays at an angle 9 1 with respect to the horizontal edge 
(X-axis, see Fig. 2). For 9i < n/2 — 9 C , conversion to a P ray can occur only at the 
vertical edge (Y-axis, see Fig. 2). The P-ray at an angle 9 2 (measured from the X-axis), can 
however reconvert to an S ray at both the vertical and the horizontal edges. If it reconverts 
at the vertical edge, the resultant reflected S-ray is again at an angle 9 X while if reconversion 
occurs at the horizontal edge, the 5*-ray (which we denote by S') is at an angle 9 3 > 9 X (see 
Fig. 1). The first possibility requires no further analysis while for the S-vay at an angle 
#3, reconversion can occur only at the horizontal edge and hence generates a P ray at angle 
9 2 . An identical process occurs if 9i > 9 C except that edges are now interchanged. For 
7r/2 — 9 C < 9i < 9 C (see Fig. 1) and for 9 X = or 7r/2, there is no conversion and the orbit 
continues with the same intensity. This exhausts all possibilities so that there are only three 
directions {9 X , 9 2 , 9 3 } that all rays originating from a given parent ray can explore. 

As k, increases, the two branches in the map of Fig. 1 move apart and the range of angles 
at which conversion can occur at either the vertical or horizontal edge shrinks. However, the 
scenario described above continues to hold whenever conversion does occur and the daughter 
rays produced explore only three momenta directions. 

As k decreases however, the branches in the map (of Fig. 1) come closer and meet at 
k c = y/2. As k decreases below k c , there is a range of angles in which conversion occurs 
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both at the horizontal and vertical edge as shown in Fig. (3). Consequently, the number 
of directions accessible to all daughter rays from a single parent increases as well. As k 
decreases further, the overlap increases and at k — > 1 + , conversion can occur at both the 
horizontal and vertical edge for any initial angle. Thus rays split with every reflection and 
the number of daughter rays grows as 2 k where k denotes the number of reflections. The 
collection of all daughter rays now explore an increasing number of momentum directions in 
sharp contrast to the case when k > k c . 

Note that in all cases, the orbits are marginally unstable as in the case when no conversion 
occurs. This can be verified by linearising the neighbourhood and looking at the eigenvalues 
of the Jacobian matrix. Periodic rays thus occur in families and their extent is limited by 
the vertices where adjacent parallel rays meet a horizontal and a vertical edge respectively 
and convert differently. 

III. PERIODIC ORBITS 

We now turn to a study of periodic orbits in such a system and provide an algorithm for 
determining them when k > k c . 

We first note that unlike a normal square billiard, each set of winding numbers (M, N) 
can correspond to more than one periodic solution. The total length, l s , of all S segments 
(at angle 6\) in a periodic trajectory is such that Is cos(6 l i) = miL, where mi is an integer 
and L is the size of the square. Similarly, the total length, ls> of S' segments in a periodic 
trajectory is such that ls> sin(^ 3 ) = n\L where n\ is again an integer. It follows then that 
the total length, l P of all P segments in a trajectory with winding numbers (M, N) is such 
that : 

l p cos(# 2 ) = [(2M - mi)L - ni L cot(0 3 )] (8) 

On equating the total projected length of the trajectory along the vertical edge (Is sin(9i) + 
l P sm(6 2 ) + Is 1 sin(0 3 )) to 2NL we obtain : 
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mi tan(0i) + [(2M - mi) - n x cot(0 3 )] tan(0 2 ) = (2JV - m) (9) 

where m\ can vary from to 2M while n\ can vary from to 2N . Note that # 2 and 8 3 can 
be expressed in terms of 8\ so that for a given value of (mi, ni), the root of eq. (Q) (if any) 
can be determined. However not all (mi, n\) admit real solution though in general any set of 
winding numbers, (M,N) allows more than one periodic orbit. The case m x = 2M,n 1 = 
corresponds to a normal periodic orbit that suffers no conversion. 

The total length of all S ray segments in any trajectory with winding numbers (M, N) 
and labelled by (mi, m) is Is + Is' — Lm\/ cos(#i) + Ln\/ sm{8 3 ) while the length of P ray 
segments put together is lp = (2M — mi — n\ cot(6> 3 ))L/ cos(# 2 ). The time period of the 
trajectory is thus (lg + ls')/ c s + lp/cp. 

Clearly, these S, S' and P segments can be arranged in several ways so that an orbit with 
a given (M,N,m,i,ni) will in general be degenerate. In order to compute the degeneracy, 
it is necessary to construct symbol itineraries that are not related by cyclic permutation. It 
turns out that four distinct symbols, S, S', P, P' together with the winding numbers (M, N) 
specify a trajectory. Here S(P) and S'(P') refer to a single stretch of S(P) ray between two 
consecutive vertical and horizontal edges respectively (see Fig. 2) though it must be noted 
that P and P' have the same angle 82 where as S and S' are at angles 8\ and #3 respectively. 
Also, the total length of all P segments is in general not equal to m 2 L/ cos(6 , 2 ) +n 2 £/ sin (#2) 
where m 2 and n 2 are integers that count the number of P and P' segments respectively. 
This is due to the fact that there exist small segments of P ray joining adjacent vertical and 
horizontal edges. 

Consider for example the case with (mi, m, m 2 , n 2 ) = (1,1,1,1) [|1^]. There are six distinct 
itineraries : (S, S', P, P'), {S, S', P', P), (S, P, S', P'){S, P, P', S'), {S, P', S', P), (S, P', P, S') 
and some of these correspond to periodic trajectories with winding numbers, (M, N) = 
(2, 2) UTB| . It is implicit however that small segments of P rays joining adjacent vertical 
and horizontal edges exist, whenever the following two symbols occur consecutively in an 
itinerary : (S, S'),(S',P),(S,P') or (P,P')- Fig. 2 illustrates this for (S',P,S,P') with 
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(M,N) = (2,2). 

Finally it is important to remark that the above symbols do not specify a trajectory 
uniquely and it is possible that more than one itinerary corresponds to the same orbit. This 
occurs when the total length of a single P segment in a trajectory is such that it hits at 
least two vertical and also two horizontal edges. 

For k < k c , the analysis gets increasingly complicated though as k — > 1 + , all daughter 
rays born out of a parent periodic orbit at an angle 9\ = tan -1 (iV/M) are eventually periodic 
with the same length (for cs = cp = 1, the time period equals the length). The degeneracy 
is thus 2 2 ( M+7V ) so that the density of orbit lengths can be expressed as 

<*(0 = ££<*(*" 2LVM 2 + TP) e 2( - M+N ^ 2 (10) 

M N 

The average proliferation rate can thus be obtained by integrating over M and N as : 

d av (l) = JdM JdN 5(1 - 2LVM 2 + A 2 ) e 2{M+N ^ 2 (11) 

= f^de J dr^ 5(1 _ r ) e ^lB2cos(e- 7 r/4)/i (12) 



where 2LM = r cos(#) and 2LN = r sin(9). The r integration is trivial and the 9 integration 
can be evaluated asymptotically for large I using the Laplace method fl4| . We finally obtain : 



= TT^lhl e " (13) 
where the growth exponent h = -^ln 2. For L = 7r/4, h = 1.248. 

As k —> oo, no conversion is allowed and the average density of periodic orbits is then 

27i7 

dav(l) = ^ (14) 

With this background, we now present some numerical results on the proliferation rate 
of periodic orbits for k > k c . 



IV. NUMERICAL RESULTS 



We first consider the 1.429 and choose cs = 1 and cp = k. Using the procedure 

outlined above, we have generated all periodic orbits that have time periods less that 7.5. 
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We then construct the staircase function, N(T) = J2i ®(T — 7$) which counts the number 
of periodic orbits with Tj < T. Fig. 4 shows a plot of In N(T) as a function of T together 
with the best fitting straight line. The fit is good indicating an exponential proliferation of 
periodic orbits. The slope, which is a measure of the growth rate, equals 1.509. 

Note that a direct comparison with the growth rate at k — > 1 + (obtained in section pTT[ ) 
is not possible due to the fact that the velocities at k — 1.429 are necessarily different. 
Assuming that all conversions are allowed, the mean velocity required to achieve a growth 
rate of 1.5 is about 1.2. 

We next consider the cases when k equals 8 and 30 and choose cs = 1/k and Cp — 1. 
The range of angles which permit conversion at k = 8 is now much smaller and a large 
fraction of periodic orbits do not undergo any conversion for the time periods considered. 
The proliferation rate is thus sub-exponential (for a plot, see Fig. 5). 

For k = 30, conversion is further inhibited as the range of angles (in which conversion is 
allowed) becomes smaller still. This is evident in Fig. 5 where we plot In N(T) as a function 
of In T for k = 8 and 30. The curve for k = 30 has a good linear fit with slope 2.22 indicating 
a power law behaviour (note that the exponent is expected to be 2 in the limit k — > oo). 
The curve for k = 8 follows this till In T ~ 2 and then increases as more converted orbits 
are included. 

This substantiates our analysis and shows that with a decrease in the average number 
of conversions as k increases, the proliferation law for short periods shows a crossover from 
exponential to power law behaviour. However for large but finite k, exponential proliferation 
is eventually expected to dominate for large T though these orbits are not easily accessible 
to computations. 

V. DISCUSSIONS AND CONCLUSIONS 

We have, in the preceding sections, analysed a simple ray-splitting system and shown that 
as the parameter k increases from 1 + , mode conversion is progressively inhibited and does 
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not take place with every reflection at the boundary. This shows up in the proliferation law 
of periodic orbits as well. Surprisingly, exponential proliferation persists for k ~ \pl where 
the number of directions accessible to the daughter rays is only three. For short periods, 
there is a crossover to sub-exponential proliferation with increasing k since conversion can 
occur in a progressively shrinking range of angles that are (nearly) parallel to the two edges 
and converted periodic orbits are thus longer on an average. 

In terms of conservation laws, the decay in intensity that accompanies splitting must 
be compensated by a faster proliferation of periodic orbits, though not necessarily of the 
exponential kind. The sum rule for ray-splitting however needs to be derived in order 
to be more specific, though heuristically, its form should be similar to Eq. (0) with each 
term having an additional factor representing the intensity loss. This is an important area 
to explore for sum rules can be put to practical use for example in checking whether all 
periodic orbits up to a certain length have been determined. 

Our computations were limited by the fact that a set of four symbols together with the 
winding numbers were necessary to label periodic orbits. We have, in each case included all 
periodic orbits with symbols strings of length 10 and the period T up to which all periodic 
orbits are available is decided by the shortest orbit with symbol string of length greater than 
10. For boundary shapes leading to hyperbolicity, the complexities increase making longer 
orbits practically inaccessible to computations. Thus, even though our results for larger 
values of k are limited to short orbits, they are significant in this light. 
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FIGURES 
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FIG. 1. The map shows the relationship between the angles of the S and P ray on conversion 
at k = 1.429. All angles are measured with respect to the X axis; i.e. horizontal edge. Notational 
details can be found in the text. 
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FIG. 2. Two unfolded periodic orbits belonging to the same family with winding number (2,2) 
and symbol itinerary (S 1 , P, S, P'). Note the small segments of P ray between S' and P as well as 
S and P'. 
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FIG. 3. The map shows the conversions possible at k = 1.05. The angles 9s and 6p denote 
respectively the angles of the S and P ray measured with respect to the X axis (horizontal edge). 
Note that there is a range of angles where conversion can occur at both the horizontal and vertical 
edges. This is referred to as an overlap. 
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FIG. 4. A plot of In N(T) as a function of T together with the best fitting straight line. Here 
N(T) counts the number of periodic orbits that suffer conversion and have periods less than or 
equal to T. The growth exponent as derived from the slope is 1.509. 
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FIG. 5. The proliferation law for k = 8 (upper curve) and k = 30 (lower curve). Note that we 
plot In N(T) as a function of In T. Also shown is the best fitting straight line for k = 30. 
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